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Spawning Models for the CPHD Filter 

Daniel S. Bryant, Emmanuel D. Delande, Steven Gehly, Jeremie Houssineau, Daniel E. Clark, Brandon A. Jones 


Abstract —In its classical form, the Cardinalized Probability 
Hypothesis Density (CPHD) filter does not model the appearance 
of new targets through spawning, yet there are applications for 
which spawning models more appropriately account for newborn 
objects when compared to spontaneous birth models. In this 
paper, we propose a principled derivation of the CPHD filter with 
spawning from the Finite Set Statistics framework. A Gaussian 
Mixture implementation of the CPHD filter with spawning is 
then presented, illustrated with three applicable spawning models 
on a simulated scenario involving two parent targets spawning a 
total of five objects. Results show that filter implementations with 
spawn models provide more accurate results when compared to 
a birth model implementation. 

Index Terms —Multi-object Filtering, CPHD Filter, Point Pro¬ 
cesses, Random Finite Sets, Bayesian Estimation, Target Track¬ 
ing, Target Spawning 


I. Introduction 

T he goal of the multi-object estimation problem is to 
jointly estimate - usually in the presence of clutter, data 
association uncertainty, and missed detections - the time- 
varying number and individual states of targets evolving in a 
surveillance scene. Commonly known detection and tracking 
algorithms for the multi-object problem include Joint Proba¬ 
bilistic Data Association (JPDA) [1] and Multiple Hypothesis 
Tracking (MHT) [2]. Relatively new is the multi-object filter¬ 
ing framework known as Finite Set Statistics (FISST) [3], [4], 
based on a representation of the target population as a Random 
Finite Set (RFS), a specific case of the more general concept 
of point process. 

Within the FISST framework, the multi-target Bayes filter 
proposes an optimal solution to the multi-object estimation 
problem; it is, however, impractical in realistic applications 
due to its combinatorial complexity [3]. Several approxima¬ 
tions of the multi-target Bayes filter have been proposed to 
circumvent this intractability, including the Probability Hy¬ 
pothesis Density (PHD) [5] and the Cardinalized Probability 
Hypothesis Density (CPHD) [6] filters. The PHD filter prop¬ 
agates the first-order factorial moment density, or intensity, 
of the multi-target RFS, representing the whole population of 
targets within the surveillance scene [5]. While inexpensive, 
the PHD filter exhibits a high variability in the estimated 
target number [4]. The CPHD filter [6] addresses this issue by 
estimating the cardinality distribution of the multi-target RFS 
in addition to its intensity. Unlike for the PHD filter, the initial 
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presentation of the CPHD filter does not include a model for 
target spawning. Target spawning refers to instances where 
a parent target generates one or more daughter targets and 
where the daughter(s) usually remain(s) in close proximity to 
the parent for some amount of time following their appearance, 
e.g., a fighter jet that launches a missile. 

Though the CPHD filter’s model for birth targets has the 
potential to address spawning targets [4], there may be cases 
where specific spawning models are more applicable. In the 
context of the tracking of Resident Space Objects (RSOs), nat¬ 
ural and artificial Earth orbiting satellites consisting of active 
spacecraft, decommissioned payloads, and debris, consider for 
example the deployment of CubeSats from a launch vehicle 
[7], [8] or fragmentation events caused by the unintentional 
[9] or intentional [10] collision of objects. Without spawning, 
the best option may be the use of diffuse birth regions, 
however, the volume of space to be filled requires a potentially 
intractable number of birth regions [11]. To improve the 
CPHD filter’s performance for space-object tracking, [12] 
presented a measurement-based birth model that leverages an 
astrodynamics approach to track initialization for RSOs. While 
such an approach may be effective for tracking spawned RSOs, 
a multi-target filter that correctly models the birth process for 
a given target is expected to provide better accuracy and faster 
confirmation of new objects. The models proposed in this 
paper allow for the development of CPHD implementations 
used for RSO tracking applications with spawning. 

The incorporation of spawning models in the context of 
CPHD filtering has previously been explored in [13], relying 
on an intuitive construction of the filtering equations related 
to the spawning models considered (Bernoulli or Poisson 
process) through a non-standard derivation procedure. In this 
paper, we propose expressions for the CPHD filter enhanced 
with various target spawning models through a standard 
derivation procedure within the FISST framework specific to 
the considered spawning model (Bernoulli, Poisson, or zero- 
inflated Poisson process). To the best of our understanding, the 
derivation of the spawning terms in [13] relies on additional 
approximations and the approach does not lead to the same 
results as those presented here. 

The structure of this paper is as follows. Section II presents 
the relevant background on point processes and functional dif¬ 
ferentiation, followed by key definitions and properties perti¬ 
nent to our results. Section III provides a detailed construction 
of the CPHD Alter with target spawning, considering several 
models of spawning processes. Section IV demonstrates the 
proposed concepts through simulation example, and closing 
remarks are given in Section V. The proofs of the results in 
Section III are given in the Appendix. 
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II. Background 

In this section, we introduce the necessary background 
on point processes (Section II-A), on Probability Generating 
Functional (p.g.fl.)s (Section II-B), on functional differen¬ 
tiation (Section II-C), and on a few properties from the 
application of differentiation in the context of point processes 
(Section II-D). 

A. Point processes 

A point process on some space X is a random variable 
whose number of elements and element states, belonging to 
X, are random. In the context of multi-target tracking the 
population of targets is represented by a point process $, on a 
single-target state space X C K'^, whose elements describe 
individual target states. A realization of $ is a vector of 
points (fi = {xi,... ,xn) depicting a specihc multi-target 
conhguration, where G X describes the d-component state 
of an individual target (position, velocity, etc.). 

A point process $ is characterized by its probability 
distribution on the measurable space (A’,B.%’), where 
'^ = Un>oX" is the point process state space, i.e., the space 
of all the hnite vectors of points in X, and B;f is the Borel 
(T-algebra on X [14]. The probability distribution of a point 
process is dehned as a symmetric function, so that the order 
of points in a realization is irrelevant for statistical purposes 
- for example, realizations {xi,X 2 ) and {x 2 ,xi) are equally 
probable. In addition, if the probability distribution is such that 
the realizations are vectors of points that are pairwise distinct 
almost surely, then the point process is called simple. For the 
rest of the paper, all the point processes are assumed simple'. 

The probability distribution P<j, is characterized by its 
projection measures P^\ for any n > 0. The n*'‘-order pro- 
jection measure P^ , for any n > 1, is dehned on the Borel 
(T-algebra of X" and gives the probability for the point process 
to be composed of n points, and the probability distribution 
of these points. By extension, is the probability for the 
point process to be empty. For any n > 0, denotes the 

^th 

-order Janossy measure [16, p. 124], and is dehned as 

4”^(Pi X ... X P„) = ^ P4^(Pai X ... X (la) 

a{n) 

= n!P4^(Pi X ... X P„), (lb) 


ways to describe the point process $. However, a measure- 
theoretical formulation provides a more general framework 
that is required to construct certain statistical properties on 
point processes that can be exploited for practical applications; 
a recent example is given in [17] for the construction of the 
regional statistics. For the sake of generality, the rest of the 
paper thus uses a measure-based description. 

Assuming that / is a non-negative measurable function on 
X, then the integral of / w.r.t. to the measure P$ can be 
written in the following ways: 


PM) 


[ f{‘f)PMp) 
Jx 

(2a) 

/ f{p)pM)dp 

J X 

(2b) 


Xn)Pt\dixi,...,Xn)) (2c) 

X / 

... ,x„)da::i ... 

(2d) 

n>0 

•. ,a:„)4"V(2;ii ■ • ■ 

(2e) 

n>0 ^ ^ 

,Xn)j^\xi, . . .,Xn)dxi . ..dXn. 

(2f) 


Throughout this article the exploitation of the Janossy mea¬ 
sures will be preferred, for they are convenient tools in the 
context of functional differentiation (see Section II-C). For 
the sake of simplicity, domains of integration will be omitted 
when they refer to the full target state space X. 

The Janossy measures can also be used directly to exploit 
meaningful information on the point process $. For example, 
central to this article is the extraction of the cardinality 
distribution of the point process, that describes the number 
of elements in the realizations of $ (see Section III): 


Example 1 (Cardinality distribution). Consider the function 
fn defined as 


fn{p) 


1, \'p\=n, 

0, otherwise, 


( 3 ) 


where Bi is in Bx, the Borel tr-algebra of X, 1 < i < n, and 
where ^{n) denotes the set of all permutations [ai ,..., (t„) 
of (1,... ,n). 

The probability density (respectively (resp.) the n*'*- 
order projection density the -order Janossy density 

fnl 

) is the Radon-Nikodym derivative of the probability dis¬ 
tribution Bj, (resp. the n*'*-order projection measure P^\ the 
n*'‘-order Janossy measure 4"^) with respect to (w.r.t.) some 
reference measure. All these quantities provide equivalent 

*An alternative construction of simple point processes as random objects 
whose realizations are sets of points ip = {xi,..., xjsf}, in which the 
elements are per construction imordered, is also available in the literature 
[6], [15]. In this context, a point process is called a RFS. 


where \ip\ denotes the size of the vector ip. The integral of fn 
w.r.t. to P(j, yields the probability p^{n) that a realization p 
of the point process $ has size n and we have, using Eq. (2) 
(see [18, p.28]): 

pM)=PMn) (4a) 

= [ Pi"'\d{xi,...,Xn)) (4b) 

= [ J^\d(^Xi,...,Xn)). (4c) 

n\ Jx" 

The function is called the cardinality distribution of the 
point process 4>. Note that the -order projection measure 
P^^ (resp. the rf^-order Janossy measure ) is not a 
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probability measure, in the general case, for its integral over 
X" yields p,i,(n} (resp. n\p^{n)). 


B. Probability generating functionals 

The p.g.fl. provides a useful characterization for point 
process theory [19] and is defined as follows. 


Definition 1 (Probability generating functional [16]). The 
probability generating functional of a point process $ 
on X can be written for any test function h G 14 (X) as^ 


G^{h)= f P<5(d(/5) (5a) 

xGip 

= ^ f H^l) ■ ■ ■HXn)J^\d{xi, . . . ,Xn)). 


(5b) 


The p.g.fl. G$ fully characterizes the point process <i>, and is 
a very convenient tool for the extraction of statistical informa¬ 
tion on $ through functional differentiation (see Section II-C). 
From Eq. (5) we can immediately write 

G$(0) = J® (=P®), (6) 

G$(l) = 1. (7) 


Operations on point processes (e.g., superposition of two 
populations) can be translated into operations on their corre¬ 
sponding p.g.fl.s. In the context of multi-target tracking, p.g.fl.s 
provide a convenient description of the compound population 
(targets or measurements) resulting from an operation on 
elementary populations. 

The superposition operation for point processes describes 
the union of two populations $i, $2 into a compound popu¬ 
lation $1 U $2, during which the information about the origin 
population of each individual is lost. 

Proposition 1 (Superposition of independent processes [19]). 
Let $1 and $2 be two independent point processes defined 
on the same space, with respective p.g.fl.s G$^ and G^j. The 
p.g.fl. of the superposition process $1 U $2 is given by the 
product 

G$ju<i> 2(^) = G$i(/i)G$2(/i). (8) 


The Galton-Watson recursion for point processes [19], [20] 
describes the evolution of each individual x from a parent 
population $p into a population of daughter individuals, 
independently of the other parent individuals but following 
a common evolution model described by a process <l’e- The 
resulting daughter population (I’d is then the superposition of 
all the populations of daughter individuals. 

Proposition 2 (The Galton-Watson recursion [20]). Let G$p 
be the p.g.fl. of a parent process <l>p on X, and let G$,,(-|x) 
be the conditional p.g.fl. of an evolution process $6. defined 
for every x € X. The p.g.fl. of the daughter process $d is 
given by the composition 

G^Ah) = ■ (9) 


C. Functional differentiation 

To make use of functionals in the derivations presented in 
Section III, we require the notion of differentials on functional 
spaces. We adopt a restricted form of the Gateaux differential, 
known as the chain differential [21], so that a general chain 
rule can be determined [22], [23]. Following this, we describe 
the general higher-order chain rule. 

Definition 2 (Chain differential [21]). Under the conditions 
detailed in [21 ], the function f on some set X has a chain 
differential Sf(x; rj) at x G X in the direction rj if, for any 
sequence r]n ^ P & X, and any sequence of real numbers 
On 0, it holds that 

df{x;p)= lim if{x + O^Pn) - f{x)) ■ (10) 

n^oo tfn 

The -order chain differential can be defined recursively 
as 


=S {S^ ^f {x-pi,...,Pn-l);Vn) ■ 

( 11 ) 

Applying n*^-order chain differentials on composite func¬ 
tions can be an extremely laborious process since it involves 
determining the result for each choice of function and proving 
the result by induction. For ordinary derivatives, the general 
higher-order chain rule is normally attributed to Faa di Bruno 

[24] . The following result generalizes Faa di Bruno’s formula 
to chain differentials and allows for a systematic derivation of 
composite functions (see [22] for an example of exploitation 
in the context of Bayesian estimation). 

Proposition 3 (General higher-order chain rule, from [23], 

[25] ). Under the differentiability and continuity conditions 

detailed in [25], the -order variation of composition fog 
in the sequence of directions at point x is given by 


5-{fog){x-,{p,)U) 


7rEn„ 


( 12 ) 


where n„ = n({l, ..., n}) represents the set of partitions of 
the index set {1,..., n}, and |7r| denotes the cardinality of the 
set TT. 

Example 2 (General higher-order chain rule). 

O 9)ix-,pi,P2) 

= S‘^f{9ix);6g{x;pi),6g{x;p2)) + Sf (g(x); pi, 772)) ■ 

^={{1}.{2}} ’^={{1.2}} 

(13) 

Applying n*^-order chain differentials on a product of 
functions follows a more straightforward approach, similar to 
Feibniz’ rule for ordinary derivatives. 


^W(X) is the space of bounded measurable functions u on X satisfying 

||u||oo < 1. 


Proposition 4 (Feibniz’ rule, from [25]). Under the differen¬ 
tiability conditions detailed in [25], the -order variation of 
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the product f ■ g in the sequence of directions (?7i)r=i point 
X is given by 

^"(/• 5 )(x;(^.)r=i) 

7rC{l,...,n} 

where = {1,..., n} \ tt denotes the complement of tt in 
{ 1 , 

Example 3 (Leibniz’ rule). 


III. The CPHD filter with spawning 

This section covers the derivation of the filtering equations 
for the CPHD filter for various target spawning processes. 
Section III-A provides a brief description of the general 
multi-target Bayes filter [3], and the principled approxima¬ 
tion leading to the construction of the original CPHD filter 
[ 6 ]. Section III-B then presents the various models of point 
processes that will be necessary for the construction of the 
CPHD filter with spawning in Section III-C. 


■ 9){x\m,V2) = 

T]i,r]2)g{x) -f 6f{x\ r]i)Sgix; 92) 

7r={l,2} •"•={1} 

+ Sf{x; g2)5g{x; 771 ) -f f{x)Sg{x; pi, 772) ■ (15) 

7r={2} •7r={0} 


D. Probability generating functionals and differentiation 
Key properties of a point process can be recovered from 
the functional differentiation of its p.g.fl.. Taking the fc‘^-order 
variation of G$(/i) in the directions pi,..., p^, we have (see, 
for example [26, p. 21 ]), 

5'=G$(/i;pi,...,Pfc) = 

7 -TTj / Y[v^ix^) WKxi) 4"^(d(a;i,...,a:„)). 

(16) 

It is then useful to consider the cases when we set h = 1 or 
h = 0, i.e. 



S'"G,s>{0',Vi, - ■ ■ ,Vk) 

= f m{xi)---Vk{xk)J^^\d{xi,...,Xk)), (17) 

7x'“ 

^'^'G<i,(l;pi,... ,pfc) 

= [ m{xi) ■ ■ ■Vk{xk)M^\d{xi,... ,Xk)), 

7x'“ 

where is the /c‘^-order factorial moment measure, de¬ 

fined as in [14, p. 111]. 

Assuming that one wishes to evaluate the Janossy and facto¬ 
rial moment measures in some measurable subsets Bi £ Bx, 
1 < i < k, then they can be recovered from Eqs (17), (18) 
by setting the directions to be indicator functions^ p^ = 1^^, 
1 < i < k, so that 

= X ... X Bk), (19) 

S^G^ (/i; 1 b, ,..., 1b JI (Si x ... x (20) 

The propagation of the first-order factorial moment measure 
- also called the intensity measure - of the multi¬ 
target point process $, in a Bayesian context, is a key 
component of the construction of both the PHD filter [5] and 
the CPHD filter [6]. The density of the intensity measure is 
called the Probability Hypothesis Density [5]. 


^For a measurable subset B G Bx, the indicator function Is is defined 
as the function on X such that la (a;) = 1 if a: G S, la (a;) = 0 otherwise. 


A. Multi-object filtering and CPHD filter 

The multi-target Bayes filter [3] is the natural extension of 
the usual single-target Bayesian paradigm to the multi-target 
case, within the FISST framework. The multi-target Bayes 
recursion at time step k consists of the time prediction and 
data update steps given as follows: 


-Ffe|fe-i(d7?|^i:fc-i) = [fk\k-i{p\p)Pk-i{dip\Zi,k-i), (21) 

Jx 


Pk{dnp\Zi,k) = 


gk{Zk\p)Pk\k-i{'^p\Zi,k-i) 

JX 9k {Zk \p')Pk\k—l {dp\Z\-k—\) 


( 22 ) 


where Pk\k-i (resp. Pk) is the probability distribution of 
the predicted multi-target process 4)fe|fe_i (resp. the posterior 
multi-target process 4>fc), Zi, 1 < i < k, is the set of mea¬ 
surements collected at time step i, Zi-i denotes the sequence 
Zi,...,Zi, fk\k-i is the multi-target transition kernel, and 
gk is the multi-target likelihood function. The multi-target 
transition kernel fk\k-i describes the time evolution of the 
population of targets since time step k — 1 and encapsulates 
the underlying models of target birth, motion, spawning, and 
death. The multi-target likelihood gk describes the sensor 
observation process and encapsulates the underlying models 
of target detection, target-generated measurements, and false 
alarms. 

The multi-target Bayes recursion is used to propagate the 
posterior distribution Pk{-\Zi.k) that describes the current 
target population based on all the measurements Zi,..., Zk 
collected so far. The CPHD Bayes recursion aims at simpli¬ 
fying the multi-target Bayes recursion by approximating the 
predicted and posterior multi-target processes as independent 
and identically distributed (i.i.d.) processes"^, a class of point 
processes fully characterized by their cardinality distribution 
and their first-order moment measure [ 6 ]. The CPHD 
filter thus focuses on the propagation of the posterior cardinal¬ 
ity distribution pk and the posterior first-order moment mea¬ 
sure Pk, rather than the full posterior probability distribution 

Pk. 

The original construction of the CPHD filter [ 6 ] does not 
consider a target spawning mechanism, and the key contribu¬ 
tion of this paper is to propose the integration of several target 
spawning models in the CPHD time prediction equation (see 
Section III-C). Note that the data update step does not involve 
the target spawning mechanism and is therefore left out of the 
scope of this paper. A detailed description of the data update 
step can be found in [27]. 


^The definition of an i.i.d. process is given in Section III-B. 
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B. Point process models 

1) Bernoulli process: A Bernoulli process $ is character¬ 
ized by a parameter 0 < p < 1 and a spatial distribution s. 
It describes the situation where 1) either there is no object in 
the scene, or 2 ) there is a single object in the scene, with state 
distributed according to s. Its projection measures are given 
by 

{ l-p, n = 0, 

ps(Bi), n = 1, (23) 

0 , otherwise. 

Proposition 5 (p.g.fl. of a Bernoulli process [4]). The p.g.fl. 
of a Bernoulli process $ with parameter p and spatial distri¬ 
bution s is given by 

G:s>{h) = 1 — pp j ft.(a;)s(da;). (24) 

2) Poisson process: A Poisson process $ is characterized 
by a rate A > 0 and a spatial distribution s. It describes 
a population whose size follows a Poisson distribution and 
whose individual states are i.i.d. according to s. Its projection 
measures are given by 

. . A” ” 

...xB^) = e-^—Y\s{B,). (25) 

n! 

Proposition 6 (p.g.fl. of a Poisson process [4]). The p.g.fl. of 
a Poisson process 4> with rate A and spatial distribution s is 
given by 


G<i>{h) = exp 


A 


h{x)s{dx) — 1 


(26) 


3) Zero-inflated Poisson process: A zero-inflated Poisson 
process 4) (from [28]) is characterized by a parameter 0 < 
p < 1, a rate A > 0, and a spatial distribution s. It describes a 
population that is 1 ) either empty, or 2 ) non-empty, with size 
following a Poisson distribution and whose individual states 
are i.i.d. according to s. Its projection measures are given by 


1 —p-bpe 


n = 0, 


Pi"\Bi X ...xBn) = { 






Z =1 


otherwise. 


(27) 

Note that a Poisson process is a special case of a zero-inflated 
Poisson process in which the parameter p is set to one. 


Proposition 7 (p.g.fl. of a zero-inflated Poisson process). The 
p.g.fl. of a zero-inflated Poisson process 4> with parameter p, 
rate A, and spatial distribution s is given by 


{h) = 1 — p -b p exp 


A h{x)s{dx) — 1 


(28) 


4) I.i.d. process: An i.i.d. process 4> is characterized by 
a cardinality distribution p and a spatial distribution s. It 
describes a population whose size is distributed according to 
p, and whose individual states are i.i.d. according to s. Its 
Janossy measures are given by 

n 

4”^ (Bi X ... X S„) = nlp{n)l[s{B,). (29) 

i=l 


Note that a Poisson process is a special case of i.i.d. process 
in which the cardinality distribution p is Poisson. 


C. Prediction step 

In this section, we propose an alternative expression of the 
original CPHD time prediction step [ 6 ] in which newborn 
targets originate from a spawning mechanism rather than 
spontaneous birth. Note that the assumptions on the posterior 
multi-target process from the previous time step, the target 
survival mechanism, and the target evolution mechanism are 
identical to the original assumptions in [ 6 ]. 

Theorem 1 (CPHD with spawning: prediction step). Assuming 
that, at step k: 

• The posterior multi-target process is an i.i.d. 

process with intensity measure pk-i, with cardinality 
distribution pk-i, and spatial distribution Sk-i, 

• A target in state x at time k — 1 survived to time k with 
probability ps^kix), 

• A surviving target in state x at time k — 1 evolved since 
time k — 1 according to a Markov transition /s,fc(-|a:), 

• There was no spontaneous target birth since time k — I, 
m Newborn targets were spawned from prior targets (see 

next page), 

then the intensity measure pk\k-i tind cardinality distribution 
Pk\k-i of the predicted multi-target process ^k\k-i ore given 
by 


Pk\k-l{-) = J [Ps,k{x)fs,k{-\x) -b Ph^k{-\x)] Pk-l{dx), 


(30) 


Pfc|fc —1 (^) — ^ ^ {bl, ... ,bn) 

i=i 


E 

m>j 


. n!(m — j)! 




, (31) 


where Bnj is the partial Bell polynomial [29] given by 


kdn,k{x\ : X 2, '' ' , X^i ) 


E 


nl 


ki-\-2k2-\ - \-nkn—'B' 

ki-\-k2-\ - \-kn—k 


cG ^ 2 ^ v. 


(32) 


and where the intensity measure ph,k ond the coefficients bi 
are the parameters of the spawning process, dependent on the 
modeling choices. Denoting ps^k{') = l~Ps,fe(’) ondph,k{') = 
1 — Ph,ki'), the parameters are as follows: 

a) Bernoulli process, with parameter ph,k and spatial dis¬ 
tribution Sb,fe-' 


Ph,k{-\x) = Ph,k{x)sh,k{-\x) 


(33) 
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and 

' j Ps,k{x)ph,k{,x)sk-i{dx), 1 = 0 , 

fj / [Ps,fc(x)pb,fe(a;)+ps.fe(a;)pb,fc(a:)] Sfe_i(dx), i = 1, 

2 Jps,kix)ph,kix)sk-i{dx), i = 2, 

,0, i>2. 

(34) 

b) Poisson process, with rate Ab,fe and spatial distribution 

Sh,kt 

Ph,k{-\x) = Ab,fc(a;)sb,fc(-|a;), (35) 

and, 

b^ = J 

X [ps,kix)Xh,kix) + ips,kix)\ Sk-iidx) 

i > 0, (36) 

c) zero-inflated Poisson process, with parameter ph,k, rate 

Xb,k, nnd spatial distribution Sb,fc-' 

Ph,k{-\x) = Ph,kix)Xh,kix)sb,ki-\x), (37) 

and, 

f Ps,k(x) [pb,k{x) + Pb,k{.x)e~^'^’>‘'^^'i] Sfc_i(da::), 

i = 0, 

/ [Ps,k 

, _ I +Ps,kix) [pb.kix) Pb,kix)e~^'^’’‘''^'^] ]sfc-i(da::), 

— \ . ^ 

I = 1 , 

X [ps,kix)Xb,kix) -\-ips,kix)]sk-iidx), 

i > 2. 

(38) 

The proof is given in the Appendix. Note that the structure 
of the predicted cardinality (31) allows for its efficient com¬ 
putation through an algorithm dedicated to the computation of 
partial Bell polynomials (see [29], [30] for examples). 


1) Filtering assumptions: We follow the usual assumptions 
of the GM-CPHD filter [27] regarding the transition process 
from time fc — 1 to time k, namely, that the probability 
of survival ps.fe is uniform over the state space X and the 
transition fs^k follows a linear Gaussian dynamical model: 

Ps,ki-) = Ps,k, (39) 

fs,k\k-li-\^) =-^i- '^PkX,Qk), (40) 

where Af{- ; m, P) denotes a Gaussian distribution with mean 
m and covariance P, Fk is a state transition matrix, and Qk 
is a process noise covariance matrix. 

Regardless of the chosen spawning model (see Theorem 1), 
we further assume that the spatial distribution of each spawned 
object Sb,fc can be described as the Gaussian mixture 

Sb.fc(-lx) = ^ (41) 

( 7 ) ( 7 ) 

where d], ^ is a deviation vector, ^ is a spawning transition 

' (j) ’ 

matrix, and Q], ^ is a spawning noise covariance matrix, for 
1 < j < Jb,k, and '^hk ~ 1 - Also, we assume that the 

model parameters pb,k, Xb,k, when applicable, are uniform 
over the state space X: 

Ph,k{') = Pb,k, (42) 

Xb,k{-) = Ab,fe. 

2) Predicted intensity: The construction of the predicted 
intensity Pk\k-i in Eq- (30) follows a similar structure as for 
the usual GM-CPHD filter [31]. Assume that the posterior 
intensity pk-i can be written as a Gaussian mixture of the 
form 

Jk-1 

Pk-i{-) = ( 43 ) 

i=i 

(j) 

where rn^_^ (resp. Pk_i) is the posterior mean (resp. covari¬ 
ance) of the j-th component of the mixture. Then the predicted 
intensity Pkjk-i can also be written as a Gaussian mixture of 
the form 


IV. Simulation 

In this section we illustrate the CPHD filter with spawning 
models through a simulation-based scenario. The Gaussian 
Mixture (GM) implementation of the CPHD filter is briefly 
described in Section IV-A, followed by a description of the 
metrics exploited for the analysis of the Alter results in 
Section IV-B. The scenario and the selection of the Alter 
parameters are detailed in Section IV-C, and the results are 
discussed in Section IV-D. 

A. The GM-CPHD filter with spawning 

Since the incorporation of spawning in the CPHD Altering 
process does not affect the data update step, we shall focus 
in this section on the speciflcs of the prediction step for the 
GM-CPHD Alter with spawning. A description of the usual 
GM-CPHD, including the implementation of the spontaneous 
birth term, is given in [27]. 


Pk\k-li') — Ps,k\k-li’) + lJ-h,k\k-li')y (44) 

where the surviving component fJ,s,k\k-i is the Gaussian 
mixture 


ps,k\k — l (■) 

with 


i=i 


; m 


U) 

s,k\k—l 


pU) 


), (45) 


^ii\k-i = ( 46 ) 

PlMk-i = Qk + FkPt^iFk, (47) 

for 1 < J < Jfc-i, and the spawning component p,b,k\k-i is 
the Gaussian mixture 


3 = k 
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with 


m 


{hi) 

b,A;|A;— 


1 


pU^i) 

-^b,fe|fe-l 


+ h') 


Vb,fe 


_l_ p(*) p(i) 


(49) 

(50) 


for I < j < Jfe-i, 1 < * < Jb.fc, and the scalar ah,k depends 
on the spawning model: 


Pb.fc, 

Bernoulli process, 


^b.fc, 

Poisson process. 

(51) 

Ph^kXh^ki 

zero-inflated Poisson process. 



3) Predicted cardinality distribution: Due to the assump¬ 
tions presented in Section IV-Al, the coefficients of the Bell 
polynomial in Eq. (31) have the simpler form 


a) Bernoulli process: 



(1 -Ps,k) (1 -Ph,k) , 

i = 0, 

Ps,k (1 - Ph,k) 4- (1 - Ps,k)Ph.k, 

i = 1, 

‘^Ps,kPh,k^ 

1 = 2 , 

0 , 

i > 2 . 


b) Poisson process: 

= Cfe [(1 - + ^Ps,k ], * > 0. (53) 

c) zero-inflated Poisson process: 




(1 Ps,k') (l Ph,k 4" Ph^kt^ ^ 

(1 - Ps,k)Pb.ke~^^-'‘Xh,k 

+Ps,k (l - Ph,k + Ph,ke~^^’>’) , 
PbjkX^j^ 6 [(1 Ps,fc) 4“ tPs,/c] ; 


i = 0, 

i = 1, 
i > 2. 
(54) 


The predicted cardinality distribution is then computed by the 
appropriate substitution of Eqs. (52)-(54) into Eq. (31). 


B. Evaluation metrics 

To compare the multi-target state representing the true 
targets in the scene - the “ground truth” - and a collectior 
of targets extracted from the filter’s output, we exploit the 
Optimal Sub-Pattern Assignment (OSPA) metric [32] for as¬ 
sessing the accuracy of multi-object filters. Given two sets 
X = {xi,... ,Xm}, Xi € X, 1 < i < m, and Y = 
{yi, ■ ■ ■ ,yn}, yj G X, l < j < n, the second-order OSPA 
distance d^\x,Y) between X and Y is defined as 

di^\x,Y) = 

'0, m = n = 0, 

rf / m V 1 1/2 

< -immy^d^‘'\xi,y^U)y+c^{n-m)] ,m<n, 

*- ^ 1—1 ^ -■ 

d!" 2 ^ (F, X), otherwise, 

(55) 

with 

d^‘'Hxi,yj) = min(c, ||xi - yj\\), (56) 


where c is the cutoff parameter, and || • || is the usual norm 
on X. The OSPA distance is such that 0 < d^\x,Y) < c; 
dh\x,Y) = 0 indicates that X and Y are identical, while 
Sf^{X,Y) increases with the discrepancies between X and 
Y, taking into account mismatches in number of elements and 
element states. 

In order to compare the true number of targets in the 
scene and a estimated cardinality distribution extracted from 
the filter’s output, we exploit the Hellinger distance [33]. 
Given two finite cardinality distributions P = {pi,... ,pk) 
and Q = {qi,... ,qk), the Hellinger distance dniP, Q) is 


dniPi Q) = 


y2\ 




(57) 


2=1 


Note that in (57), the coefficient 1 /a/ 2 is included in order 
to scale the Hellinger distance such that it is bounded as 
0 < dniPjQ) < 1; duiPjQ) = 0 indicates that P and Q 
are equivalent, where as dniP^Q) P and Q become 

increasingly dissimilar. 


C. Scenario and filter setup 

A point [x,y,x,y] of the single-target state space X C 
describes the position and velocity coordinates of an object 
in a square surveillance region of size 2000 m x 2000 m. The 
simulated multi-target tracking scenario consists of one scan 
per second for 100 s, and up to seven targets evolving in the 
region with constant velocity. Two targets are present at the 
beginning of the scenario and each spawns targets at different 
times: target 1 spawns two additional targets at f = 15 s 
and target 2 spawns three additional targets at t = 25 s. All 
spawned targets have a lifespan of 60s. Pig. 1 shows the 
trajectories of the targets cumulated over time, while Pig. 2 
illustrates these trajectories and the collected measurements 



0 500 1000 1500 2000 

X[m] 


Pig. 1: Target trajectories. A circle “O” indicates where a 
trajectory begins, and a square indicates where a trajectory 
ends. The large square indicates the limits of the sensor’s PoV 
and the large dashed circle represents the 90% confidence 
region of the Gaussian component of the spontaneous birth 
model. 























(a) x-axis 



(b) ^-axis 


Fig. 2: Collected measurements (gray crosses) and target 
positions (black lines). 


The probability of survival ps^k (39) is constant throughout 
the scenario, and set to Ps,fe = 0.99. The target motion model 
/s,fc|fc-i (40) is set as follows: 


Fk = 


I 2 

O 2 


AI2 

I 2 


Qk — ^ 1 / 


^ 1 

-^i-2 




(58) 


where A = Is, a,y = 5ms“^, and 1„ (resp. 0„) denotes the 
n X n identity (resp. zero) matrix. 

The sensor’s probability of detection is uniform over the 
sensor’s FoV, and set at a constant value of 0.95 throughout 
the scenario. Each target-generated measurement consists of 
the target’s coordinate position with an independent Gaussian 
white noise on each component, with a standard deviation 
of 10 m. Spurious measurements are modeled as a Poisson 
point process with uniform spatial distribution over the state 
space and an average number of clutter per unit volume of 
12.5 X 10“®m“^, that is, an average of 50 clutter returns per 
scan over the surveillance region. 

For the sake of comparison, the usual GM-CPHD f 
[27] with spontaneous birth and no spawning is impleme 
as well. The spontaneous birth model is Poisson, wit 
constant rate of 0.025 per time step (which yields, over 
100 s of the scenario, an average of 2.5 newborn targets 
each parent target). The spatial distribution is modeled wi 
single Gaussian component, centered on the sensor’s Fo^ 
illustrated in Fig. 1. 

The spatial distribution of the spawning (41) is identical for 
the three considered models. We assume no spawned target 
deviation vectors, and a standard deviation of 12 units is set 
on each component of the spawning noise covariance, i.e. 


Fh,k = 


I 2 

O 2 


O2 

I 2 


dh,k — 0; Qh,k — 


CTbl 2 


O 2 


(59) 


where 0 denotes the null vector in X, 
12ms-i. 


O 2 CTbl2 
CTb = 12 m, and = 


The parameters of the three spawning models are set as 
follows. The zero-inflated Poisson model assumes one spawn¬ 
ing per parent target during the scenario with an average of 
2.5 daughter targets per spawning event, thus pb./c and Ab are 
set to 0.01 and 2.5, respectively. Relative to the zero-inflated 
Poisson model, the Poisson model is set to yield a similar 
spawning intensity thus its Ab,/c is set to 0.025, whereas the 
Bernoulli model is set to yield a similar spawning frequency 
so its pb.fc is set to 0.01. These parameters are also presented 
in Table I. 


TABLE I: Spawn model parameters. 

Model pb.fc Ab.fe Mb,fe(-k) 

Bernoulli 0.01 - 0.01A(' ; x, Qh,k) 

Poisson - 0.025 0.025Af(' ',x,Qh,k) 

zero-inflated Poisson 0.01 2.5 0.025Af(' ;3;,Qb,fe) 

It is interesting to note that neither the Poisson nor the 
Bernoulli models are equipped to capture the nature of the 
spawning events occurring in this scenario, since, per construc¬ 
tion, the Poisson model is a poor match for spawning events 
occuring at unknown dates and the Bernoulli model is a poor 
match for spawning events creating more than one daughter 
target. The zero-inflated Poisson model possesses a greater 
flexibility and should be able to cope with a wider range of 
spawning situations; in any case, it is expected to yield better 
performances on the scenario presented in this paper. 

To maintain tractability, GM components are truncated 
with threshold T = 10“®, pruned with maximum number 
of components Jmax = 100, and merged with threshold 
[/ = 4 (see [31] for more details on the pruning and merging 
mechanisms). Additionally, the maximum number of targets 
is set to Amax = 20 to circumvent issues with infinitely tailed 
cardinality distributions [27]. 

D. Simulation results 

The proposed spawning models and the birth model are 
implemented with the GM-CPHD Alter, and compared over 
500 Monte Carlo (MC) runs of the multi-target scenario 
decribed in Section IV-C. 



Pig. 3: MAP estimate of the number of targets (averaged on 
500 runs). 


The MAP estimate of the number of targets is plotted in 
Pig. 3, along with the true number of targets in the scene. 
The results suggest that the spawning models provide a better 
estimate of the number of targets and, in particular, converge 
faster to the true number of targets following the appearance 
of new targets in the scene. This is expected, because the 
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scenario does not feature any spontaneous but only spawning 
related births, and thus in this context spawning models are ; 
better match than the birth model. 

Among the three spawning models, the zero-inflated Poissoi 
converges the fastest following the appearance of new targets 
while the Bernoulli model converges the slowest. This i 
expected, for the zero-inflated Poisson model provides thi 
best match to the spawning events occurring in this scenario 
Note in particular that the Bernoulli model may not consider 
the appearance of more than one daughter per spawning 
event, and must therefore stage the multiple-target appearance 
across several successive time steps; in other words, thi 
Bernoulli is ill-adapted to “busy” events where targets appea 
simultaneously. Note also the slight overestimation shown b; 
the Poisson model when the true number of target is stable 
Per construction, the Poisson model is well-equipped for thi 
simultaneous appearance of an arbitrary number of spawne( 
targets at any time step, but it fails at coping with “quiet' 
periods where no spawning occurs because, unlike the zero- 
inflated Poisson model, it does not temper the Poisson-driven 
spawning with a probability of spawning. In other words, the 
Poisson model is ill-adapted to the spawning events shown in 
this scenario. 

Note that all models - spawning and birth - follow the 
same mechanism for target deaths and yield much closer 
performances when target disappearances occur. 




Fig. 4: OSPA distance (averaged on 500 runs). 


Similar conclusions can be drawn from the comparisons of 
the OSPA distances shown in Fig. 4. All models show error 
spikes at times of spawning (t = 15 s, f = 25 s) and death 
(f = 76 s, f = 86 s), however, the spawning models recover 
more quickly than the birth model, and have consistently lower 



(a) Predicted cardinality 



(b) Updated cardinality 


Fig. 5; Hellinger distances (averaged on 500 runs). 

The quality of the estimation of the number of targets 
proposed by the four models is further illustrated in Fig. 5, 
where the Hellinger distance between the cardinality distri¬ 
bution propagated by each model and the “ideal” cardinality 
distribution (i.e., a distribution in which all the mass is 
concentrated on the true number of targets). 

The results in Fig. 5 allow a more refined analysis of 
the proposed models. All the models yield poor estimates 
immediately after a change in the true number of targets but 
the zero-inflated Poisson model converges the fastest following 
a target birth/death and it converges to the best estimate during 
periods where the number of target is stable. The Poisson 
model converges faster than the Bernoulli model, but to a 
worse estimate; this is expected, since the Poisson model is 
ill-adapted to “quiet” periods while the Bernoulli model is ill- 
adapted to “busy” events (see discussion above on Fig. 3). 

As expected, the updated cardinality distributions are con¬ 
sistently more accurate than the predicted cardinality distribu¬ 
tions since they benefit from the processing of an additional 
measurement batch. 

V. Conclusion 

The motivation for the work presented in this paper is the 
resolution of multi-object detection and tracking problems in 
which newborn objects are spawned from preexisting ones. 
To this end, the construction of a CPHD filter in which the 
appearance of newborn targets is modeled with a spawning 
mechanism rather than spontaneous birth is proposed, based on 
a principled derivation procedure within the FISST framework. 

A GM implementation of the CPHD filter with spawning 
is then presented, considering three different models for the 
spawning mechanism based on a Bernoulli, a Poisson, or a 
zero-inflated Poisson process. The three resulting filters are 

^Recall from Eq. (57) that the Hellinger distance dn is such that 0 < 
dn < 1 . 


errors. 
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then illustrated, analyzed, and compared to a usual CPHD filter 
with spontaneous birth but no spawning, on the same simulated 
scenario involving two parent targets spawning a total of 
five daughter targets. Results show that a spawning model, 
appropriately chosen for a given application, can provide better 
estimates than a spontaneous birth model. 
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Appendix 


A. Proof of theorem 1 

Eor the sake of simplicity, the time subscripts will be 
omitted throughout the proof when there is no ambiguity. Also, 
we will denote by ps (resp. ^d) the function 1 — ps.fc (resp. 
1 -Pd.fc)- 

1) Predicted p.g.fl.: Let us focus first on the p.g.fl. Gfc|fc-i 
of the predicted multi-target point process $fe|fe_i. Each parent 
target in the population, represented by the prior point process 
<i)fe_i, generates daughter targets in the predicted population 
in two ways: 

« a daughter target stemming from the (eventual) survival of 
the parent target, represented by a survival point process 

« a population of daughter spawned from the parent target, 
represented by a spawning point process <i)b. 

Using Eq. (8), and denoting by Gg (resp. Gb) the p.g.fl. of the 
survival (resp. spawning) point process, we can describe the 
evolution of a parent target with state a; G X with a compound 
process with p.g.fl. 

Gc{h\x) = Gs{h\x)Gb{h\x), (60) 

and exploiting the Gabon-Watson equation (9), we may finally 
write 


Gkik-iW — Gk-i{Gc{h\-)) (61a) 

= Gfc_i(Gg(/i|-)Gb(/i|-)). (61b) 

2 ) Predicted intensity: Let us now focus on the expression 
of the predicted intensity iik\k-i- For that, let us fix an 
arbitrary measurable subset B G Bx. The expression of 
the intensity evaluated in B can be recovered from the first 
derivative of the p.g.fl. Gk\k-i using Eq. (20): 

Pk\k-i{B) = ^Gfc|fc_i(/i; (62a) 

= b(Gfe_i(Ge(/l|-));lB)L=l (62b) 


Using the definition of the p.g.fl. (5a) then yields 

Pk\k-i{B) = s( [ Gc(/i|x) Pk-i{d(f);lB] 

\ ^ rV L rp^fZir, -I / 


= / nGc(/i|x);li 

^ \ ^cz,r, 


Pk-i{d(p)] 1 b 

Pk-iidif) (62d) 


h=l 

(62c) 


xG(p 


h^l 


Erom the product rule (14) it follows that 

Pk\k-l{B) = 


LL 


SGc{h\x; 1 b) 




n^c(iix) 

h—l xG‘f:> J^-| 
x^x 


Pk-i{d(p) 


= 1 


(62e) 
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Using the product rale (14) on Gc(-|x) = Gs(-|a;)Gb(-|a;) then Exploiting Eq. (28) yields 
yields 




5Gs{h\x] Gb(l|a;) 




+ Gs(l|x) (5Gb(/i|a;; 

= 1 


Pk-i{dip) 

(62f) 


Using Eq. (20) we introduce the intensity ps (resp. pb) of the 
survival (resp. spawning) process and we obtain: 


Mb(S|-) 

= JGb {h\-;lB)\k=i 
= <5(pb(-)+Pb(-)exp 



(65a) 


/i(x)sb(dx|-) —1 



/ h=l 

(65b) 


= Pb(-)^b(-)'^^y /i(x)sb(da;|-) - 1; 1 b 

Ab(-)(y Sb(da:|-) - 1 


X exp 


=0 


= Pb(-)Ab(-)sb(-B|-). 


(65c) 

(65d) 


(J-klk-liB) 



[Hs{B\x) + fj.h{B\x)] Pk-i{d(p) (62g) 


Which becomes, using Campbell’s theorem [34, p. 271]: 

Mfe|fe-i(S) = y [fJ.s{B\x) + fj.h{B\x)] /ifc_i(dx). (62h) 

Note that the validity of the expression of the predicted 
intensity above is not restricted to specific models for the 
prior process <i>fc_i. As such, the construction of the predicted 
intensity is identical in the case of the PHD filter with 
spawning (see Mahler’s original proof in [5]). Let us now 
focus on the explicit expression of the intensity measure fis- 
Since the survival process is assumed Bernoulli with parameter 
Ps{-) and spatial distribution /s(-|-)^ we can exploit Eq. (20) 
to retrieve the intensity ps through the expression of the p.g.fl. 
Gs given by Eq. (24): 


3) Predicted cardinality: Let us now focus on the expres¬ 
sion of the predicted cardinality pk\k-i- From Eq. (4) the 
cardinality distribution of an arbitrary point process can be 
retrieved through its Janossy measures; let us then compute the 
predicted n*^-order Janossy measure evaluated at the 

neighborhood of a collection of n arbitrary points yi,... ,yn- 
Using Eq. (19) yields 

= S^Gkik-i{h;Uy,,...,UyJ\^^^ (66a) 

= S (Gfc_l(Gc((l|-)); Idyi, ■ • ■ , ldy„)|^^Q (66b) 

Applying the general chain rule (12) then gives 
= (Gc(/i|-);(<5'“'Gc(/i|-;(ld,JiG^ 

(66c) 



= -Ps(-)+ft(-) y h{x)fs{dx\-)-,lB 


= P.{-)fs{B\-). 


(63a) 


h=l 

(63b) 

(63c) 


Let us now focus on the explicit expression of the inten¬ 
sity measure /ib of the spawning process, depending on the 
modeling choices. 

a) Bernoulli process with parameter Pb(-) and spatial dis¬ 
tribution Sb(’l’)- 

Using the same construction as in Eq. (63) we have immedi¬ 
ately 


Developing the predicted p.g.fl. Gk-i through Janossy mea¬ 
sures with Eq. (2) then gives 




7rennm>|7rP ' 


h—0 


X J7 Gc(0|xi)J^™](d(a:i,... ,Xm)). 

i=\Tr\ + l 

(67) 


Since the prior process is assumed i.i.d., we can substitute the 
expression given by Eq. (29) to the prior Janossy densities 
Jk^i and obtain 


dh{B\-) = Ph{-)sh{B\-). (64) 


b) zero-inflated Poisson process with parameter Ph{-), rate 
Ab(-) and spatial distribution Sb(-|-)- 


Ji'ik-Myi’-'-’yn)) 


= E E 

ttEIIti m>|7r| 


ml 


(m- |7r|)! 


p{m)G^{d{yi,...,yn)), ( 68 ) 
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where 


Therefore, Eq. (72) simplifies as follows 


C'tt( d(2/i, • ■ •, 2 /ti)) 

r r 


h—0 


C^{d{yi,.. ■,yn)) = Gs{0\x)Gb{0\x)s{dx) 








/i =0 


Now, from the derivation shown in Eq. (63), we see that; 

l-Ps{x), z/= 0 , 

5''''''Gs{h\x] (ldyjiei')|;j . = >( ps{x)fs{dyr\x), v = {i}, 


0 , 


Therefore, Eq. (70) simplifies as follows; 


<5l"lGb(/z|a;;(ld,;,),6.) 


h^O 


X Gc(0|xi)]^s(da;i) (69a) 

2^|7r| + l i—1 

i-|7r| 

Gc(0|a;)s(da;) 


^ n ( / Gc{h\x;{ldy,)i(^u>)\t,^os{dx)] (69b) 

/..C-TT V / 


Recall from Eq. (60) that Gc{h\x) = Gs(/i|a;:)Gb(/i|a;); using 
the product rule (14) on Eq. (69b) then yields 


/ r V m-|7r| 

G^{d{yi,... ,yn)) = f / Ps(a;)pb(x)s(da:) j 
X n ( / ft(x)pb(x).b(dt;.|x).(dx) 

+ J Ph{x)ps{x)fs{dyi\x)s{dx) 

X n ( /Ps(x)pb(x)/.(dt/.|x)sb(d%»s(dx) 


Ps{x)ph{x)fs{dyj\x)sh{dyi\x)s{dx) 


no 


(74) 


UJ^TT 

\ui\>2 


s(da:) I. (70) 


IH>i. 

(71) 


m—|7r| 


Substituting Eq. (74) to Eq. ( 68 ), we may finally retrieve the 
scalar ptik-iin) through Eq. (4); 


Pklk-lin) 

= (15a) 

° E E |d)! '’*’"*‘’”~'"nVh (15b) 

7i-Gn„m>|7r| 5 ' '' ia;G7r 

where the coefficients bi are defined by 


G^{d{yi,...,yn)) = (^J Ps{x)GhiO\x)s{dx) 

X n ( / (ldyJiGta;)|;j^oS(da;) 

luGtt ^ 

+ Ps{x)fsidy,\x) 

Idj/^ GUJ 

X (5l‘^l-iGb(/i|a;;(ldy,)jga;\{i})|;,^o®(‘^^))- ^ 2 ) 

We shall now detail the expression of Eq. (72) depending on 
the modeling choices for the spawning process. 

a) Bernoulli process with parameter ps(’) and spatial distri¬ 
bution /s(-|-): 

We may draw similar results from the derivation shown in 
Eq. (71); 


b. = 


' f Ps(x)pb(x)s(dx), 2 = 0, 

/ lps(x)ph(x) + Ps(x)pb(x)J s(da;), i = 1, 

2 / ps(a;)pb(a;)s(da;), i = 2, 

. 0 , z>2. 


(76) 


Using the definition of the Bell polynomial (32) then yields 
the desired result. 

b) zero-inflated Poisson process with parameter Pb(’)’ rate 
Ab(-)> and spatial distribution Sb(-|-)- 
Applying the chain rule (12) to the p.g.fl. (28) yields 


(l-pb{x), i/ = 0, 

S Pb(a:)sb(dyj |a;), = {j}, 

[O, |Z 2 |> 1 . 

(73) 


(5l^lGb(/z|a;;(ld,,),g.)|^^o = 

f 1 - pb{x) + pb{x)e~^'^^^\ 


[Pb{x)e njgi.Sb(d?/j», |z^| > 0. 


’ ( 77 ) 
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Therefore, Eq. (72) simplifies as follows 

C’tt (d(2/l; • ■ • ; Un)) = 


J Psix) 

n 

{i}€7r 


Phix)+phix)e 


s(dx) 




ujGtz 

|c^|>l 


' Psix)phix)e ^*’('')Ab(a;)sb(dyi|a:)s(dx) 
Ps(a;)(pb(a:) +Pb(a:)e“^*’^“^)/s(dyi|a:)s(dx) 


^ n ( f Psix)phix)e ^'’^“Ub(a:)l“l ]^Sb(dyi|a:) 

ljGtz ^ ^ 


s(da;) 


s(da;) ). (78) 


+ J Ps(x)pb(x)e-^"(")Ab(x)l“l-i 
X ^ fsidyi\x) Sb(dyj|x) 

iG.uJ ^ j^uj 

j¥^i 

Substituting Eq. (78) to Eq. ( 68 ), we may finally retrieve the 
scalar Pfe|fc_i(n) through Eq. (4); 
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where the coefficients bi are dehned by 


(79a) 

(79b) 


/ps(x) [pb(a;)+Pb(a;)e s(dx), 

i = 0, 

/ [psix)pbix)e~^'^^^'>Xhix) 

+Psix) [phix) +pb(a:)e“^*>(^)] ]s(da;), 
bi= < ^ 

i = 1 , 

/ Phix)Xi^^ix)e-^^^^'> 

[Ps(a:)Ab(x) + ipsix)] s(dx), 

i > 2. 

(80) 

Using the dehnition of the Bell polynomial (32) then yields 
the desired result. 
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